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Abstract 

We consider the theory of the glass phase and jamming of hard spheres in the large space dimension 
limit. Building upon the exact expression for the free-energy functional obtained previously, we find 
that the Random First Order Transition (RFOT) scenario is realized here with two thermodynamic 
transitions: the usual Kauzmann point associated with entropy crisis, and a further transition at 
higher pressures in which a glassy structure of micro-states is developed within each amorphous 
state. This kind of glass-glass transition into a phase dominating the higher densities was described 
years ago by Elisabeth Gardner, and may well be a generic feature of RFOT. Micro states that are 
small excitations of an amorphous matrix - separated by low entropic or energetic barriers - thus 
emerge naturally, and modify the high pressure (or low temperature) limit of the thermodynamic 
functions. 
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I. INTRODUCTION 

The Random First Order Transition (RFOT) scenario for glasses introduced by Kirkpatrick, Thirumalai and 
Wolynes [1-4] proposes that the glass transition is represented - at least at the mean-field level - by a freezing 
transition similar to the one of the Random Energy Model [5] , described mathematically by a one-step replica sym- 
metry breaking (1RSB) ansatz [6]. Although this scenario was first based on an analogy with spin glasses [1, 2], it was 
quickly realized in a pioneering work by Kirkpatrick and Wolynes [7] that the glass transition of <i-dimcnsioiial hard 
spheres in the limit d — > oo could be described within the same framework. Much later [8, 9], similar results were 
obtained using the replica method, which also allowed for a detailed description of the glass phase and in particular 
of the jamming point where the pressure of the glass becomes infinite, corresponding to its close packing - which was 
therefore called glass close packing (GCP) and is closely related to the random close packing concept introduced by 
Bernal much earlier [10]. 

The main advantage of the replica method is that it allows for a unified treatment of both the glass and the jamming 
transitions, within a simple static RFOT scenario, and it also allows for a partial understanding of dynamical aspects. 
Moreover, there is hope that the replica method can provide an exact result in the limit d — > oo. A first step in 
this direction was performed in the first paper of this series [11], where we have shown that the thermodynamics of 
hard spheres in the limit of high dimensions may be exactly obtained from the knowledge of the distribution of the 
two-point correlation function between states, encoded in the Parisi parameter. In the same paper it was shown that 
once a 1RSB ansatz is made, one recovers exactly the Gaussian replica free energy that was used in Ref. [8, 9] to 
derive estimates of the various transitions that characterize the RFOT scenario at the 1RSB level. 

In this paper we take a second important step, by investigating the stability of the 1RSB solution towards further 
levels of replica symmetry breaking. We find that for higher pressures, well above the RFOT one, there is a second 
transition (a so-called Gardner transition [12]) leading to a somewhat different physics in that limit, and in particular 
around the jamming point. We believe that this physics is intimately connected with the peculiar mechanical properties 
of jammed states of hard spheres, that have been recently characterized in much detail [13, 14]. 

The rest of the paper is organized as follows. We start our presentation by a general discussion of the RFOT 
scenario, of its connection with the physics of jamming, and of the main new features that are due to the presence of 
the Gardner transition. This discussion is reported in Sec. II and it is for the moment mostly speculative, although 
some parts of it have been previously studied in spin glass models. Next, we present our new results, which constitute 
a first important step to substantiate this picture for hard spheres in the d — > oo limit. In Sec. Ill we provide a 
proof of the correctness of the Gaussian ansatz for a generic form of the overlap matrix, extending the main result of 
Ref. [11]; in Sec. IV we recall a few important results of Ref. [9, 11] that are directly needed here; in Sec. V we present 
our main results for the Hessian matrix of the 1RSB solution and in particular its so-called replicon eigenvalue, that 
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FIG. 1: A sketch of the phase diagram. 



is responsible for the instability of this solution (i.e. the Gardner transition); in Sec. VI we discuss the cubic terms 
in the expansion around the 1RSB solution and from them we extract the dynamical exponents that characterize the 
glass transition; in Sec. VII we present an approximate calculation to obtain an order of magnitude for the Gardner 
transition pressure in finite d; in Sec. VIII we summarize and draw our conclusions. 

II. A GENERAL RFOT SCENARIO FOR THE GLASS AND JAMMING TRANSITIONS 



As is by now well known, Kirkpatrick, Thirumalai and Wolynes' scenario for the liquid-glass transition involves 
a first point at which the equilibrium state fractures into an exponential number of ergodic components: this is 
the dynamical temperature Td (or pressure Pd, see Fig. 1), also called Mode-Coupling temperature because in low 
dimensions it can be computed using Mode-Coupling theory. The ergodic components are only truly dynamically 
separated in the mean-field limit, while in a realistic short-range finite-dimensional situation the system is still ergodic, 
although the dynamics slows down. As the temperature is lowered, or the pressure increased, the number of metastablc 
states contributing to equilibrium diminishes, until a point is reached where the equilibrium system is left with only 
the deepest amorphous states: this is the Kauzmann point, beyond which the thermodynamics stays dominated by 
(or "frozen in") those states. From a purely equilibrium point of view, one may picture the situation at P > Pk (or 
T < 2~k) as in the sketch of Fig. 2, with widely separated states of "size" q, defined for example as: 



with a, b two copies (replicas) of the system and k a vector of length comparable to the inverse of the inter-particle 
distance (alternative definitions of q are possible, see [15] for a review). A more precise way of stating the same thing 
is to introduce the effective potential V{q) [16], which counts the logarithm of the number of configurations having 
correlation exactly q with a "reference" equilibrium configuration. One obtains a picture as in Fig. 2, where one sees 
that configurations are either close, with overlap q = <7ea (with probability 1 — to), or very far away - in other states, 
with overlap q = - with probability m, corresponding to the two minima in the effective potential and to the two 
peaks in the Parisi order parameter P(q). Note that P(q) may in general only be non-zero where V(q) takes the 
minimal value, and that <7ea plays the role of the Edwards- Anderson order parameter [6]. 

This construction concerns equilibrium configurations, but may be generalized to describe metastable states [17, 18] 
by choosing the "reference" configuration, rather than from equilibrium, from a system perturbed by a small "pinning 
field", itself thermalized at a higher temperature T' = T/m. Technically speaking, following Monasson [18], this 



A. The generic phase diagram of RFOT models 
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FIG. 2: A sketch for equilibrium Tg < T < or Pq > P > Pk- Top: a cartoon of the free energy and its minima of width 
qea- Bottom: the Parisi order parameter P(q) and the effective potential V{q). 
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FIG. 3: Effective potential associated to higher metastable states. Metastability disappears at the threshold level. 

amounts to the following calculation: one considers m weakly coupled replicas at temperature T, takes an equilibrium 
configurations of one of the replicas as the reference configuration, and couples to it an additional replica who is 
forced to stay at distance q from it. Then one computes the free energy of the additional replica, and averages it 
over the other m. This amounts to fixing the Parisi parameter m, rather than choosing the value that maximizes the 
free energy. In this way, one obtains a bistable form for the effective potential, up to a threshold value T[ h = T/mth 
at which the minimum close to <7ea disappears (Fig. 3), and at precisely the threshold level, the stability matrix 
corresponding to the minimum at <7ea develops zero modes, signaling the fact that the states close to the threshold 
level are marginal. This shows up within the replica scenario as the vanishing of the "replicon" eigenvalue [19], and 
within the Thouless- Anderson-Palmer approach [6] as the free-energy Hessian developing zero eigenvalues. A crucial 
result of Rcf. [20] is that the out of equilibrium aging dynamics happens exactly at this threshold level, and it exploits 
the marginality of the threshold states to explore phase space. 

Let us now consider higher pressures, or lower temperatures. In most systems, there is a second transition discovered 
by Gardner [12, 21, 22] years ago, at which each state itself breaks into smaller substates. The sketch one usually 
makes is as in Fig. 4. To be more precise, we consider what happens with the effective potential and P(q) for an 
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FIG. 4: A sketch below the Gardner transition at equilibrium, T < To or P > Pq. Top: a cartoon of the free energy. Bottom: 
effective potential V(q) and Parisi distribution P(q) beyond the Gardner point. The parameter m is related to the probability 
of being in some state with q* < q < qea, which is given by 1 — m. 

equilibrium configuration beyond the Gardner point. The situation is depicted in Fig 4: there are many configurations 
at all distances between q* and qea- the state of size q* has fractured into many subcomponents of smaller sizes. 
However, going away from a configuration, up to correlations smaller than q*, one finds big barriers and no states, 
up until completely different states, having minimal overlap are reached. The Parisi parameter to is now related to 
the probability of being in some state closer than q* , i.e. within the "metabasin" [23] : this probability is given by 
1 — to. This fracturing of a state into many smaller ones also happens at the level of metastable states [21, 24]: 
there is a line in the phase diagram where all states undergo a Gardner transition (Fig. 1). Metastable states may be 
found as above [21], by considering m as a free parameter. However, in the more complex regime beyond the Gardner 
transition it is not clear how to compute the threshold level that will dominate the dynamics [21, 22]. It is possible 
that the threshold level could be identified by looking at the stability properties of the fluctuations at the level of q* , 
but this need to be clarified. See [25] for some initial steps in this direction. 

B. Vibrational modes and dynamics close to jamming 

A system of hard spheres when compressed suddenly ends up in a configuration that is blocked, with the exception 
of a small percentage of "rattlers" that are free to move within a "cage" made by their neighbors. A mechanically 
stable system of hard constituents such as this may be hypostatic, hyperstatic or isostatic, depending on whether the 
number of contacts is less, more, or precisely just what it takes to guarantee mechanical stability. Because the system 
is prepared with a rapid compression, it seems unlikely that it will be hyperstatic, because if at some time during the 
compression it reaches stability, it is unable to move on to create further, redundant contacts. The hypothesis that 
is usually made is that, forgetting the rattlers, the rest of the system is precisely isostatic, the assumption being that 
there are no "rattling clusters" other than isolated rattling particles. For a detailed discussion of the fundamental role 
of isostaticity in jammed packings see Refs. [13, 14]. An additional assumption that seems to be justified in practice is 
that isostaticity is "irreducible" , in the sense that there is no subset of particles that is separately isostatic: if in such 
a system a contact is broken, then by definition all the particles in the system eventually become mobile. Clearly, 
this is a critical situation. Indeed, it has been proposed in Ref. [13, 26, 27] that these packings are marginally stable 
from a mechanical point of view, and from this most of the anomalous scalings that are found numerically have been 
derived analytically. In particular, the criticality manifests itself in the spectrum of vibrations D{oj) at densities just 
below jamming [28, 29], which has a general shape as in Fig. 5, where one has to distinguish two features: 

• There is a branch of higher frequency modes, whose lowest frequency is a;*. The frequency ui* goes to zero as 
the pressure goes to infinity. 
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FIG. 5: A schematic picture of the spectrum close to jamming. 



• Within the gap < to < to* there are the acoustic modes, which exist even at finite pressures. Moreover, in 
Rcfs. [26, 28-32] it was shown that the softer modes do not look like plane waves, therefore acoustic modes are 
mixed with other kinds of soft modes. 

In the rest of this section we will argue that the Gardner transition provides a natural explanation for the presence 
of soft modes at u; < to* . These modes should appear at all pressures beyond the Gardner transition. However, the 
connection between the soft modes observed in Rcfs. [26, 28-32] and the ones associated to the Gardner transition is 
not clear for the moment. 

Consider the squared displacements Ai(t,t') = \xi(t) — Xi(t')\ 2 , where i labels the N particles of system, and 
its average over particles A(t,t') = N" 1 J2iLi Ai(£, t'). In the following we will assume that the system has been 
prepared by some rapid compression at time t — 0, in such a way that if t > t' > and il is large enough, the system 
is stuck into a glass state. The mean square displacement is given by the average of the squared displacement over 
the dynamical process, 



1 N 

A(M') = <A(M')> = £<M*) " ^(Ol 2 ) (2) 



and the variance of the squared displacement defines the so-called four-point susceptibility 



XA {t,t')=N 



A(i,t') 2 )-(A(M')) 2 



1 [<M0 ~ Mt')\ 2 \ Xj (t) - x 3 (t')\ 2 ) - (\ Xl (t) - Xl (t')\ 2 )(\x 3 (t) x 3 (t')\ 2 )] 



(3) 



N 



These definitions can be made more precise to take into account the presence of rattlers, we refer the reader to Ref. [33] 
for a detailed discussion. The "cage size" is the limit 

A 2 (oo) = lim lim A(t,t') (4) 

t— t' — >oo t'^-co 

where e oo' stands for times £, t' as large as the lifetime of the state. At pressure P, the natural scale of the displacements 
is 1/P, therefore it is convenient to introduce a scaled cage size as 



POO 

Aoo = P 2 A 2 (oo) oc / doj 
Jo 



(5) 



and the last relation is derived in Rcfs. [27, 33]. For d > 2, this quantity is finite for finite P, because D(to) ~ io d 1 
in the low frequency acoustic branch, however it diverges as P — > oo because to* goes to zero and the integral is 
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dominated by duj ~ D< £i (while the integral in < uj < uj* docs not contribute to the divergence). The 
fluctuations of the cage size yield the four-point susceptibility (or "spin glass susceptibility" ) [33] 

/>0O j-y/ \ 

X4 (oo)= lim lim X4(t,*') ~ / dw-M- ( 6 ) 

t — V— >oo t — too Jq UJ 

On the one hand, from the theory of the Gardner transition, we expect ^4(00) to diverge there, and to stay infinite up 
to infinite pressure. In fact, one may convince oneself that this is so just by considering the curvature of the effective 

potential above and below the Gardner transition, where d = 0. On the other hand, we may look at this from 
the point of view of normal modes: we split (6) in a contribution above, and one below uj*: 



D(u) D{uj*) 

X4 (oo) = / duj -AJ. + _A_Z (7 ) 



It has already been remarked in Ref. [33] that in three dimensions even the acoustic modes will make the integral in 
< ui < uj* diverge for finite P. However, the effect of acoustic modes shows up in Xi{^i t') om y at very long time 
differences t — t' 3> 1 so that in Ref. [33] it was shown that the regime of t — t' ~ 1/ui* gives a good definition 
of the four-point susceptibility. In our large-dimensional case (actually, for all d > 4), the density of acoustic modes 
is negligible, but we still expect that the first term in (7) diverges below the Gardner transition. The conclusion 
seems to be that there are other soft modes (below uj*) that do not contribute to the linear susceptibilities or to the 
short-time t — t' value of Xii^^'), but dominate the limit of lim t _ t /_ i>00 \4(t,t'). Although it is tempting to identify 
these modes with the ones observed in Refs. [26, 28-32], more work is needed to clarify the connection. 



C. Out of equilibrium dynamics in the Gardner phase 



The out of equilibrium dynamics of this system has not been solved, but from the structure of states one may already 
guess its main features. Below the Gardner transition line, the slow compaction (aging) dynamics should proceed 
close to the threshold level, defined as described above as the one where the stability at the level of q* is marginal. 
The relaxation process can be seen as a dynamical exploration of phase space starting from a completely correlated 
state (q = 1) down to a completely decorrelated state (q = 0). The relaxation should be fast from correlation q = 1 
down to (JeAj and then proceed - in a progressively slower way as the system ages - down to q* , and from there to 
zero. The fluctuation-dissipation properties may be studied by considering a system with hard spheres in a thermal 
bath of temperature unity, subjected to a pressure P generated by either a piston or by coupling to the radii of all 
spheres. The response and correlation functions are as described in Refs. [34, 35]: the response is computed from 
the staggered displacement R(t,t') = Yli^( x i) induced subjecting particles to random unit fields £j with an energy 
term £7/j e id = h^ i £ t i8xi, per unit of h. The conjugate correlation may be taken to be the quadratic displacement 
A(t,t') defined above. Response and correlations may be put together in a plot, which should look as in Fig. 6 for 
long times t' , the time t > t' being used to produce a parametric plot of R versus A. Here, A* and Aea are the 
values corresponding to the correlations q* and <7ea- At every time, the first barriers encountered are the small ones 
close to qEA- beyond the Gardner transition, where small states are separated by relatively small excitations because 
there are states at all distances q with q* < q < Qea- This might help explain the paradox that the path between 
these small states is mainly along the flattest vibrational modes - as found in Ref. [36] - while this is not what one 
expects for the large-scale relaxation within a supercooled liquid, at least within the RFOT scenario. 



D. Low temperature excitations 



A long standing problem in the physics of glassy and amorphous materials is the low temperature behavior of their 
specific heat and thermal conductivity, which turns out to be quite different to that observed in crystals. These features 
are all the more intriguing because they tend to be quite universal for all amorphous materials. The usual explanation 
for this phenomenon is to attribute it to localized quantum- mechanical two- level tunneling systems [37, 38]. These 
models assume that there are particles or groups of particles that evolve and tunnel in random local potentials. These 
potentials arc usually proposed phcnomcnologically, although they are of course generated by the same interactions 
that produced the amorphous solid in the first place. Furthermore, these simple localized clusters or particles will be 
coupled, and their interactions might generate collective effects. A proposal to take these features into account [39] is 
to consider a system of strongly coupled localized deformations, which one may assimilate as "spin"-like excitations, 
and to assume that they have essentially random interactions - with long range, partly because elasticity is long range, 
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FIG. 6: A sketch of an FDT plot for an aging system below the Gardner transition. The plot is composed of two straight and 
a curved segment. 

and partly to make the system solvable. One obtains in this way a "spin-glass" of deformations, with elementary 
excitations which one may calculate and which tend to be universal because of their collective nature. Note that this 
way, phenomenology has been pushed one step up, to the effective interaction of excitations. 

Quite clearly, the mechanisms that generate the coupling between low-temperature excitations, and the one respon- 
sible for the amorphous matrix on which they live, are one and the same. One would thus expect the same theory 
to explain both features. In the context of this paper, it is very tempting to interpret the large valleys (metabasins) 
of size q* as being the amorphous structure, and the excitations of all sizes between q* and (?ea as a "spin-glass" of 
small excitations within that amorphous structure. Formally this is clearly so, a fact that was already recognized by 
Gardner in her original paper, where she showed that the transition is essentially that of the Sherrington-Kirkpatrick 
model within each large state. More recently, the analogy between jammed packings and the Sherrington-Kirkpatrick 
model has also been underlined [40], and the idea that there is a deep connection between the marginality of jammed 
packings and low-temperature anomalies in glasses was also proposed by S. Nagel (e.g. in his talk at the ACS meeting, 
Philadelphia, August 2012). A possible difference may be noted with respect to Rcf. [39]: here the spin-glass transition 
need not (and in general will not) coincide with the liquid-glass transition at which the amorphous matrix is formed. 

III. REPLICATED ENTROPY OF INFINITE-DIMENSIONAL HARD SPHERES 

The above discussion provides several important motivations to look for an instability of the 1RSB solution in 
particle systems, akin to the Gardner transition of spin glasses [12]. Additional ones will be given by the more 
technical discussion that we now start, see Sec. IV B. We will show that a Gardner instability indeed happens in hard 
sphere systems in the limit d — > oo (and probably also in finite dimensions within the mean- field RFOT scenario). 

We will consider a system of N hard d-dimensional spheres with unit diameter, enclosed in a volume V, hence at 
density p = N/V. The packing fraction is <p = 2~ d pVd, with Vd = 7r d / 2 /r(l + d/2) the volume of a sphere of unit 
radius. In the first paper of this series [11], we derived an exact expression for the replicated entropy of this system 
in the limit of large dimension d — > oo. The result is obtained by first writing the entropy in a manifestly rotationally 
and translationally invariant form, and then performing a saddle point evaluation in the limit d — > oo. Within replica 
theory the resulting entropy is a function of the density function p(q), where the matrix q is a m x m symmetric matrix 
that encodes the overlaps q a b between different replicas. The main result of Ref. [11] was that a Gaussian assumption 
for p{q) gives the exact result for the entropy, i.e. for all thermodynamic properties of the system, exclusively in terms 
of q a b (with J2a=i lab = for all 6, because of translational invariance). In other words, no higher order parameters 
qabc,qabcd--- are necessary in the large dimensional limit. 

The proof of Ref. [11] was restricted to the 1RSB form of q a b- In this section, we will extend the results of Ref. [11] 
to obtain the replicated entropy as a function of the overlap matrix, without making any assumption on the RSB 
structure. We will start by deriving the Gaussian replicated entropy for a generic overlap matrix, and then show that 
this form coincides with the exact result. Obviously, in this section we will often make reference to Ref. [11], which 
we encourage the reader to consult before looking to the rest of the section. Another option is to skip this section 
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and take the result - as expressed by Eqs. (15) and (16) - for granted. This will be the starting point to study the 
stability of the 1RSB solution (and much more) in the following. 



A. Gaussian ansatz for a generic overlap matrix 



We have to parametrize a generic Gaussian form of p{q), or equivalently p(u), where q a b = u a ■ Ub and u a are the d- 
dimensional vectors corresponding to replica displacements, with u = {u\, • • • , u m }. We can choose a parametrization 
in terms of a m x m symmetric matrix A such that Y^a=i ^ ab ~ f° r a ^ Calling A m,m the (m — 1) x (m — 1) 
matrix obtained from A by removing the last line and column, the most general Gaussian form of p(u) is 

p{ u) = -e-^?-\A—rJu a .u b (8) 

(2nY m - 1 W 2 det(A m ' m ) d / 2 

which is normalized according to p = J T>up{u) and T>u = fn d S(J2 a u a )du\ ■ ■ ■ du m . The parameters A a b are interpreted 
as 

(u a ■ Ub) = - I Vup(u)u a ■ u b = dA ab , (9) 



for a, b e [l,m - 1], while (u a ■ u m ) = - Y.T=i ( u a 1 u b) = A am and (u m ■ u m ) = Y^a™ ( u <> 1 M &) = A mm- 

The saddle point value of q, that dominates all the integrals over p(q), is obtained as follows. We start from the 
normalization condition (note that a complete derivation of J(q), that was not reported in Ref. [11], is reported here 
in Appendix A) 



/m / m \ 
^ n s [Yi qab ) ^ 

a=l \b=l ) 



dqJ{q)p{q) oc / dq I \5 I > ' q ab ] et( d - m ) logdot 9 m ' m -^r~ 1 (^ m ' m L^ ^ ( 10 ) 



and maximizing the exponent for d — > oo leads, for a, € [l,m— 1], to (q m ' m )ab — (■^ m ' m )ab nence 9a? = dA a b, 
consistently with Eq. (9). 

To compute the replicated entropy using the general Gaussian ansatz we start from Eq. (45) of Ref. [11], which 
gives the following expression for the replicated entropy: 

S[p(q)]/N = l-\ogp(q s n-2 d - 1 ^T^2q^ , (11) 

where the function T is given in Eq. (37) of Ref. [11]. The ideal gas term is 

1 - log p(q sp ) = 1 - log p + dlog m + (m ~ 1)d + {m 2 1)d log(27r) + i log dct(i m - m ) . (12) 



The interaction term is 



Therefore 



M = 1 - \ogp+ rflogm + {m 2 l)d \og{2ne) + ± log dct(i m <™) - 2 d ~\T (j^2A 

W 2 A 

In order to obtain a simple limit d — > oo, it is convenient to define a matrix d = jyiA and a reduced packing 
fraction (p = 2 d ip/d. With this choice we have 

s[&] = ^ = 1 - log p + dlog m + {m 2 1)d log(2neD 2 /d 2 ) + ^ logdct(d m ' m ) - ~(pT (2d) . (15) 



The matrix d is a variational parameter and is therefore determined by maximization of the entropy. Defining 
F'aiffl) = ifo i the equation for d is 

(^"T 1 = 2^(2a) . (16) 
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Eqs. (15) and (16) provide the expression of the Gaussian replicated entropy and will be the starting point of all our 
calculations. 

Note that the entropy of the equilibrium glass is obtained by optimizing s[a]/m, given in Eq. (15), with respect to 
the matrix a and of m. Let us call a* and to* the optimal values, a* being the solution of Eq. (16). The reduced 
pressure p = j3P/ p of the equilibrium glass is given by 



Pglass(f) = - 



(p ds[a*] 
to* dip 



1 

m* 



1 + -$F{2&*) 



(17) 



This result shows that the pressure diverges whenever to* — > as p ~ 1/m*. Hence, the density at which to* — > 
defines the jamming point [9]. 



B. Exact computation 

We now show that Eqs. (15) and (16) can be equivalcntly obtained by an exact evaluation of the saddle point 
equations derived in Eqs. (64) and (65) of Rcf. [11]. In fact, we can make use of Eqs. (65) and (39) of Ref. [11] to 
obtain a closed self-consistent equation for q sp , which as before is the point where the argument of the integral 

~ « m / m \ 

p= / dqJ(q)p(q) = e\n d C m , d dq\{8 K> h e # logdet (lg) 
^ ^ a=l \6=1 / 

is maximum (sublcading terms for d — > oo have been neglected) . Taking the derivative with respect to q and computing 
the result in q = q sp leads to the equation 

D \ r ; m , m) -l = 2 ^ a L d. fp \ _ (1Q) 



Clearly, defining a — jjzq sp this equation is equivalent to Eq. (16). 
Evaluation of Eq. (18) at the saddle point gives the equation for A 

p = e x m d C m , d e |logdet«— -^(2^«-) ; (20) 

from which, using Eq. (78) of Ref. [11], we obtain 

-A = - logp + d log to +\{m- 1) log + d - logdctr p;m ' m - dtpT ( 2 JpV SP ) ■ ( 21 ) 

Combining Eqs. (64) and (65) of Ref. [11], using Eq. (21) and recalling the definition a = -jyzq sp we have 



( ^ 2fP 



d 



= 1 - logp + dlogm + -(to - 1) log ^— j + - logdet q sp ' m ' m ~ [^ P ) ( 22 ) 

= 1 - logp + dlogm + ^(m - 1) log (j^-^j + \ logdet a m ' m - ^yF (2a) 

which coincides with Eq. (15). This completes the proof of the exactness of the Gaussian ansatz for the computation 
of the entropy. Note that as already observed in Ref. [11] this does not imply that the Gaussian form (8) can be 
used to compute correlation functions (that encode structural properties), because the equivalence is only correct at 
the saddle point level for the entropy: this is consistent with the numerical observation of a non-Gaussian cage shape 
obtained in Ref. [41]. A computation of the cage shape is in progress and will be hopefully reported in future papers 
of this scries. 
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IV. 1RSB SOLUTION 



A. Reminder of the 1RSB solution 



in Rcf. [8, 9, 11] wc studied the 1-step replica symmetry breaking (1RSB) ansatz, which amounts in this formalism 
to assuming that all replicas arc equivalent [18]. For completeness let us recall here this result, which correponds to 
the simple choice 



1R.SB A i x 
*ab = A I °ab - 



(23) 



— ,2 

with A = jjiA. Within this ansatz, A is the "cage radius", as it is proportional to the long time limit of the mean 
square displacement in the glass [9]. Note that we use a small hat for matrices, while the wide hat just denotes 
reduced scalar variables. Hence A is a matrix while A is a scalar, and they should not be confused. 

Using the relation log det{L4((5 a f, — l/m)] m > TO } = (m — 1) log A — logm and the results of Sec. VIIB of Rcf. [11], it 
is easy to check that Eq. (15) reduces to the result of Ref. [11] for the 1RSB entropy, which is 



d, (m-l)d (m-l)cL / 2irD 2 A \ „ - n 

S [a 1RSB ] = 1 - logp + - logm + I 2 ; + I 2 J log — ^— - -tp[l - g m {A)\ , 



G m {A) = 1 —m 



dX 



/2tt 



=e 2 



(24) 



1 / / V2A-X 

5 I + erf — 



The equation for A is derived by optimizing the above results, leading to 



m — 1 
A 



ip- 



s d 9 m {A) 
dA 







1 



a dg m {A) 



1 



dA 



FmiA) 



(25) 



The function J-" m (j4) introduced here should not be confused with the function !F{a) introduced before. These are 
different functions, the first acts on a scalar while the second on a matrix. 

The physical consequences of this expression for the entropy have been derived in Ref. [8, 9], where the expressions 
of the dynamical transition density, the Kauzmann transition density, and the GCP density have been derived, with 
the scalings sketched in Fig. 1. Furthermore, at the level of the 1RSB solution, we know that A* ~ m* so we conclude 
that the cage radius vanishes as A* ~ 1/p. As a consequence of this scaling, the scaled cage size Aoo introduced in 
Eq. (5) is found to diverge as ~ p 2 A* ~ p, as noted in Ref. [33]. 



B. Inconsistencies of the 1RSB solution 



The 1RSB predictions for physical quantities were carefully compared with numerical results, around both the 
glass and the jamming transitions [9, 41-44]. Despite the good overall agreement with numerical data, one expects, 
as described above, that a Gardner transition to a full replica symmetry breaking scheme is generic. Furthermore, 
several inconsistencies have been found close to the jamming transition, at very high pressure: 

• The 1RSB solution predicts the existence of jammed packings with density ipj in the interval 2 -dd ~ @i ^ 
[6.26, log d] [9]. However, only the packings with (pj ~ logd are isostatic, with each particle in contact, on 
average, with z — 2d other particles. Instead, the packings with (pj of order 1 arc found to be hyperstatic with 
z > 2d, which is, as mentioned above, unexpected and inconsistent with numerical results. 

• In the glass phase, the exact relation between the pressure p and the contact value y(ip) of the pair correlation, 
p = 1 + 2 d ~ 1 tpy((p) [45], is violated. In particular, it is found that when cp — > ipj, p ~ d(fij/((pj — (p), consistently 
with numerical results, while 

VM = 2'-ifo,-- V ,) X l-2i-'dM ' (26) 

where the first term is the one that is consistent with the scaling of the pressure. Hence, the correct relation 
between p and y(ip) is recovered only if 2 1 ~ d d/tpj = 2/<Pj <C 1 when d — s- oo, which again suggests that the 
1RSB solution is inconsistent when (p is of order 1 and might be stable only when <p 3> 1. 
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• The scaling at large (reduced) pressure p of the cage radius A (the long time limit of the mean square displacement 
in the glass) predicted by the IRSB solution is A ~ p^ 1 , while the marginal stability argument of Ref. [26, 27] 
predicts that A ~ p~ 3 ^ 2 , which has been confirmed numerically in several studies, e.g. Refs. [26, 27, 33]. 
This exponent controls all the other exponents that characterize criticality at the jamming transition [33] and 
is directly related to the anomalous soft vibrational modes that appear at jamming [13, 14, 27, 33], hence 
reconciling the theoretical prediction with the numerical results is of extreme importance. 

• Other exponents that characterize the structure at jamming, for instance the famous (almost) square-root 
singularity in the pair correlation function [44, 46-48], are not reproduced by the IRSB solution, at least at 
the Gaussian level (a more detailed calculation of the structure functions based on the non-Gaussian theory 
developed in this series of papers is in progress and will hopefully be reported in a future paper). 

All these considerations suggest strongly that the IRSB solution is unstable, at least when pressure is large enough 
and (p is of order 1. They provide further motivations to study the stability of the IRSB solution [12], which is the 
subject of the next section. 



V. SECOND ORDER EXPANSION: THE HESSIAN MATRIX 



Here we obtain the expansion around the IRSB solution at the quadratic order. For this we need to consider a 
more general ansatz or expand around the IRSB solution, which is made possible by the general expression of the 
entropy obtained in Eqs. (15) and (16). The quadratic expansion of the entropy around the IRSB solution provides 
a stability matrix whose eigenvalues allow one to determine the stability of the solution. We will find that, as it 
happens in a Gardner transition [12, 21, 22], the IRSB solution becomes unstable when pressure is large enough. In 
finite and arbitrarily large dimensions, this happens for all tpj. However, for d — > oo, the so-called glass close packing 
(GCP) [9] which is the densest amorphous packing and has <p ~ logd becomes stable again, suggesting that the IRSB 
predictions for jamming arc still approximately useful as a starting point, but should be corrected to take into account 
its instability. 



A. General structure of the Hessian matrix 



Because the matrix a should have the sum of the elements of every column and every row equal to zero, we can say 
that the independent entries arc the elements above the diagonal of the matrix, provided that the matrix is symmetric 
and the diagonal is fixed in such a way that the constraints on the sum over the elements in a row or in column is 
satisfied. Hence, in the following we denote as 5/5a a <b the derivative taken with respect of the element a a \, with 
a < b which is assumed to be the only independent clement (hence its variation induces a variation of ctb a and of the 
diagonal elements a aa and abb)- Taking into account all this, we define the Hessian matrix as 



a^b,c^d 



2 A 2 S 2 s[a] 
d Sa a< b5a c< d 



Mi 



SacSbd + SadSbc 



M 2 



Sac + Sad + he + Sbd 



M 3 



(27) 



where the replica structure is a consequence of the structure of d 1RSB . Although the matrix M is defined by the 
above equation only for a < b and c < d, we will define it for convenience also for a > b and c > d assuming that it is 
symmetric (hence the notation M a ^b,c=id)- The prefactor 2A 2 /d is chosen for later convenience and is positive, hence 
it does not affect the sign of the three different eigenvalues of the mass matrix, which are 



X R = Mi , 

A L = Mi + (m - 1)(M 2 + mM 3 ) 
X A =M 1 + ^(m-2)M 2 . 



Defining the "entropic" and "interaction" terms 



M 



(i) 



S 2 T[v] 



M 



(£') 



a^b\c^d 



Sv a< bSv c <d 
S 2 

Sa a< b5a c< d 



(28) 



v=2a 1RSB 

logdet(a m < m ) 



(29) 
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wc have 



Mi = A 2 m[ E) - AspPM™ 



(30) 



We now compute these two terms separately. 



B. The entropic term 



To compute the entropic term it is convenient to introduce a shorthand notation ft = a m ' m . Let us also use indices 
i,j, k, ■ ■ ■ for f3 to highlight that they run from 1 to to — 1. In the 1RSB solution /3jj B = A(8ij — 1/to) has the same 
form of a but on the reduced (to — 1) x (to — 1) space. Then we have 

1 



1RSB\-1 

Jij 



A 



(Sij + 1) 



Moreover we have 



Using the standard formula 



we have 



6Jk 
Sa a 



<b 



SiaSjb + SibSj a - 5ij(5 ia + Sib) 



8 log det(/3) 



-^logdct(/3)=£/^/^=Tr 
oa a<b e -r1 l J oa a<b 



Sa a 



<b 



From the dchnition j3j3 1 = I where I is the identity matrix we have 

S 







Sa c<d 



oa c<d 



It follows that 



M 



(E) 

a<b:c<d 



Sa a< bSa c< d 
1 



log det(/3) 



Tr 



Sa c< 



Sa a<b Sa c<d 



A 2 



from which it follows that 



[2 (S ac Sbd + SadSbc) + (Sac + S a d + 5 b c + Sbd)] 



m[ e) = m^ e) ^-^- 

A 2 

M (E) = . 



(31) 



(32) 



(33) 



(34) 



(35) 



(36) 



(37) 



C. The interaction term 



We now consider the interaction term. Hence we need an expansion of the function T(v) around the 1RSB solution. 
Recall that v is a to x to symmetric matrix such that the sum of the elements in each row and in each column is 
zero. Starting from the results of Section V of Ref. [11], we can write explicitly the function !F(v), introducing 



m-dimensional vectors x a such that x a ■ Xb = v a b, as 



2tt) t 



■ cxp 



1 

-min \e + x a 

Z a 



lim 

n->0 



E 

Til , . . .,Tl m n a— n 



lim 

n->0 



E 



ni ,■■■,'"'771 ; 



ni! . . . n m ! 



m! . . .n m ! 



cxp 



cxp 



^E>*i 2 nE 



l,m 

1 v-^ n a n b 



2^ n< 

a.b 



~X a ' Xb 



m l,m 



2 " n 

a=l 



2 " 



"Uab 



i. 1RSB value of T 
First let us compute J 7 on the 1RSB solution where 



vT B = 2A(6„ 



1 



Let us call <5 a;m i n a function that is equal to one only if a is such that A a is the minimum among all the {A a }: 
other words min a A a = X a S a min- Then we have 



^[t3 1rsb ] = lim Y 

711 ,....n m ;X]L Tl n a— n 



m\ . ..n m \ 



cxp 



lim e A 

n-»0 



E 



ni! . . .n m ! 



exp 



lim e 

Tt->0 



E 



Tii! . . .ra m ! 



E 

2=1 
< 

771 

n 



-1 l: m 

1 v-^ n Q n 6 



_ \ " n °iRSB 

a— 1 a.b 



,1RSB 



~ u ab 



a=l 

dA„ 



exp 



m \ 2 m \ 



a=l 



a=l 



n 



dA a 

-V/27T 



exp 



-E^f-™ minA --^ 



/(n^)-p[-if(*. 

\a=l / L o=l 



dAi -k(M+V2A 
=e 



V27T 

dAi -^(Ai+V2A 



where 



/2tt 



erf (x) = —j= I dye y 

\/7T 



dA 
■ < 

V2^ 



-7S 



We also define for later convenience 



Note that 



dyy k e~ 



'2n Jx 



e a (x 
Qx(x 
e 2 (a- 

63(2; 

e 4 (x 



e(-x/V2) 

e-^ 2 /V27 
a;6i(x) + e(-x/V2) 

ei(x)(2 + x 2 ) 

x{3 + x 2 )e 1 (x) +36 (x) 
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and so on. Eq. (40) provides the derivation of the interaction part of Eq. (24) [9, 11]. 



2. Monomials of n 



We now want to expand the quantity J- around the 1RSB solution. The part of the Hessian matrix coming from 
the interaction term is: 



M {1) - — \v 1RSB ] - lim 

°< 6;c<d ~ 5v a<b Sv c<d [V l -n™ ^ m m ! . . . n,- ' 

ni a ...,n m ; / /' n a —n 



S 2 T 



E 



f(n a ,n b )f(n c ,n d )exp 



m 9 

nf 



where the functions / are at most quadratic: 



n a n b n 2 n\ n a n b 



f(n a , n b ) = — + — % L 



By introducing the following notation 



(O) = lim V 



nil . . . n m l 



Oexp 



ni,...,n m ;^]™ n a —n 

we have that the Hessian matrix is given by 



-^E?f 



a=l 



M 



a<b;c<d 



1 ( n 2 n c ,n 2 n d ,n\n c ,n\n d ,n a n 2 ,n a n 2 d ,n b n 2 ,n b n 2 d ^ 



,n a n c n d , ,n b n c n dy ,n c n a n by ,n d n a n b , 
-2 I - 21 I ) - 2 ° - 2 ( ° 
x n 3 ' x n 3 ' x n 3 ' x n 3 ' 

1 ( , n i n l\ , , n l n d \ + , n l n l \ + i n l n d \ _ y n l n c n d ^ _ 2 , n 2 b n c n d _ ^ 



ir' 



rr 



ntn a n b 



ir' 



2l n d n a n b + ^ n a n b n c n d , 
n 4 n 4 

Hence we want to compute averages of monomials of the {ro a }, which can be written as follows: 



fl fl I 11. — Vfl 71 & — > 



— j -. n ai ■ ■ ■ n ak cxp 

Hi,...,n m ;^J? 1 n a —n 



h n 



n->0 
k c -A 



dA a 
v/2^ 



(~l) fc e 
(2A) fc / 2 



si n 



dA„ 



E 



ni! . . . n m ! 



E™! 



2tt 



E 



■ exp 



9 s 



2AJ2 



n a \ a 



(44) 



(45) 



(46) 



(47) 



ni! . ..n m l dX ai ■ ■ -dX a , 



exp 



^ n 

a— l 



(48) 



(2A) fc /2 
1 

(2A) fe / 2 
1 



V a 

n 



9* 



dA„ 



exp 



— V 2A min A a 



dA n 



:>-i--<a = l 2 - 



9* 



<9A Ql • • • d\ ak 



,-pm A£ 



exp 



oE ( Aa + ™<5a 



>E™i 



9* 



(2j4)*/a \ 9A ai ■ • • d\ ak 

where the definition of the average has been changed to 



(o>=/ ii| oh p 4e( a «+^ 



(49) 
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The factor in parenthesis is a polynomial in {A a }, hence we now want to be able to write averages of monomials of A. 
Using the replica symmetry of the average over {A a }, we need in particular the following objects: 

B a b = ( " 9 b \ = — ^ ( — S a b + A a Afc) 

\ n 2 / 2A 

T a bc = ( a t / = — ^tt: (Sbc^a + S ac Xb + 5 a b\ c — A a AbA c ) 

\ / (2A) 3 / 2 (50) 



<n a n b n c n d \ 1_ 
) 2 

Scd^a^b — <5bdA a A c — 5 a( /AbA c + A a Af,A c A(i) 



Aobcd = / a - C d \ — ^ (^aft^ed + 0~bc.5 a d + ^ac<5 c d — ^afcA c Arf — Sbc^a^d — Sac^b^d 

\ n 4 / (2A) 2 



and the Hessian matrix is 

M a<b;c<d = 4 + B ad + B bc + B bd ) 

— j (7a ac + T aa( j + Tbbc + Tfetd + T acc + T a dd + T bc c + T b dd — 2T ac d — 2T bc d — 2T abc — 2T a bd) (51) 

~l~ T (A aacc ^aadd 



3. Monomials of A 

We therefore need to compute several monomials of the {A a }, which are listed in the following. Calculations follow 
Eq. (40) and it will be convienient to define one more average over A as 



Then wc have 



dA -hi A+V 2A 

</(A)> = -4=e 2 \ V J /(A) (52) 



2tt 



(l)^[« 1RSB ] = (me (A) m - 1 ) 
<A{) = (\ k 60(A)" 1 - 1 + (m - 1) e fe (A)e (A) m - 2 ) 
(A1A2) = (A fc e / (A)6 (A) m - 2 + A i e fc (A)6 (A) m - 2 + (m - 2)6 fe (A)e ; (A)e (A) m - 3 ) 

(A^A^Ag ) =(A fc e ; (A)e„(A)e (A) m - 3 + A i e,(A)e„(A)e (A) m - 3 + A"e fe (A)e ; (A)e (A) m - 3 

+ (m - 3)e fc (A)e ; (A)6„(A)e (A) m - 4 ) 

(Ai A2A3 a^) =(A fc e ; (A)e„(A)e p (A)e (A) m - 4 + A i e fc (A)e„(A)e p (A)e (A) m - 4 
+ A"e fc (A)e i (A)e p (A)e (A) m - 4 + a" e fc (A)e i (A)e„(A)e (A)™- 4 
+ (m- 4)e fe (A)e i (A)e„(A)e p (A)e (A)™- 5 ) 



(53) 



and so on. 



^. JTie structure of the mass matrix 
Thanks to replica symmetry the Hessian matrix has only three independent matrix elements. These are 
^i ( 2;i2 = ^(-Bii + #12) + (T ua - T111) + J(2A mi + 6A 1122 - 8A 1112 ) , 

M i ( 2;i3 = \( B n + 3B i2) - 5(2111 + T n2 - 2T 123 ) + i(A im + 3A 1122 - 4A m2 ) , ( 54 ) 
M i ( 2; 3 4 = Bu - 2(T U2 - Tus) + (Auaa - 2A U23 + Ai 234 ) . 
It is however convenient to write the matrix in this form: 

M [ a<b;c<d = MP S « cSbd + SadSbC + M 2 (7) ^c + Sad + Sbc + S^ + ^ ^ 
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FIG. 7: Phase diagram of the 1RSB solution in the (m,tp) plane (to compare with Fig. 1, recall that rougly speaking pressure 
is the inverse of m), including its instability. The black line at m = 1 corresponds to the liquid phase. The blue line is the 
dynamical line (pd{m) at which a non-trivial 1RSB solution appears. The red line is the instability line ipG(m) where the 
replicon vanishes. The red dashed line is its asymptotic behavior for small m and large (p. 



where 



Mf ) = 2A4 / : ) 12 - 4Mg 13 + 2a4 7 : ) 34 = 2A 1122 - 4A 1123 + 2A 1234 

(56) 



M 2 (/) = 4Af 1 ( £ 13 - 4M$3 4 = B n - B 12 - 2T m + 6T 112 - 4T 123 



+Aim — 4Am 2 — An 22 + 8An 23 — 4Ai 23 4 
M 3 (J) = Mg ] M = B 12 - 2T ua + 2T 123 ) + A 1122 - 2A 1123 + A 1234 

The above equation, together with Eqs. (50) and (53), give the complete expression of the interaction part of the 
Hessian matrix. 



D. The replicon 



The stability of the 1RSB solution depends crucially on the replicon eigenvalue, Xr = Mi = A 2 m[ E ^ 
Collecting all the above results and simplifying some terms we get 



±(pA 2 MY 



X R = -4-2£A m (A) 

A m (A) = (e (X) m - 5 [e>i(X) 2 - A9i(A)e (A)][(2 - 2A 2 )6 (A) 2 + (to- 4)6i(A) 2 + (6 - m)A6 (A)e 1 (A)]) 



e (Ar 



61(A) 

e (A) 



- A 



61(A) 
e (A) 



(2 - 2A 2 ) + (m - 4) 



61(A) 
6 (A) 



+ (6-m)A 



61(A) 
6 (A) 



(57) 



We can compute this numerically on the 1RSB solution, where A is the solution of Eq. (25), to find the point where 
Xr = and the 1RSB solution becomes unstable. The instability curve in the (to, (p) plane, which corresponds to the 
line ^g(to) on which the replicon vanishes, is reported in Fig. 7. 

Asymptotically we obtain (pQ(m) ~ m -1 / 2 for m — > 0. To explain this we must investigate the asymptotics of the 
function A(m, A) when both m and A are small. It is convenient to use Eq. (25) to eliminate tp instead of A. Doing 
this the equation for the instability becomes 



2F m {A) = -A m (A) , 



(58) 



which must be solved to obtain Ac(m) and then (pq{tti) using Eq. (25). We want to show that Ag(iti) ~ m 2 , and 
that this implies ipci 171 ) ~ to," 1 / 2 . 
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First of all let us examine the asymptotics of the different terms for large and positive A. For A — > oo we have 
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6i(A) 

e (A) 



(2 - 2A 2 ) + (m - 4) 



er(A) 



©o(A) 

It will be convenient for the following to define 



+ (6 - m)A 







A^ + A^ 




6i(A) 




6 50 


e (A) 


f A 1 ~ a7 + 


9x(A) 


m 


6m -4 


e (A) 


+ A^ + 



(59) 



A 6 



£(A) = 
£(A) ~ m 



9i(A)V 9i(A) 
e (A)y 6 (A) 



2 - 2A 2 + (m - 4) 



9i(A) 

e (A) 



+ (6 — m)A 



9i(A) 

e (A) 



2m 13m - 4 56 - 112m 
A 5 " + A 1 + A6 



(60) 



Now we expand Eq. (58) at small A. Let us recall that Q m {A) = 1 — m (0o(A) m 1 ). From this we obtain 



Q m (A) = G 1 (m)VA + G 2 (m)A + 

G 2 (m) 



Gx{m) 
G 2 (m) 



1 - m 2 
dA 



'2-k 
dA 



A + ... , 

1 — m 

A2 / 2 e (A) m - 1 (-AV2) 



(61) 



e - A2 / 2 e (A)"^ 1 (A 2 -l) , 



and similarly (the fact that the horrible integral corresponding to A m (A — 0) is exactly can be proven by a scries 
of integrations by parts) : 



A m (A) = Li(m)V A + L 2 (m)A + 
dA 

Li(m) 
L 2 {m) 



V2ir 
dA 



Hence Eq. (58) becomes 



with 



AA 1 (m) + AA 2 (m) 



e~ x / 2 9 (A) m - 1 £(A)(-AV2) , 
e- A2 / 2 e (A) m - 1 £(A)(A 2 - 1) . 

Ai(ra) 



(62) 



A 2 (m) 



(63) 



. , . Gitm) , . 



A 2 (m) 



1 — m 
2 G 2 (m) 
1 — m 



(64) 



+ L 2 (m) 



Asymptotically for small m, the integrals in the above expressions can have different behaviors, depending on the 
behavior of the integrand for large A when m — >• 0. In fact, if the integrand decays faster than 1/A, the integral 
is well defined and has a finite limit for m — > 0. In the opposite case, the integral is divergent and the divergence 
is dominated by the large A behavior: in this case one has to analyze the possibly divergent part to determine the 
behavior of the integral at m — > 0. In the case of Ai(m), thanks to a subtle cancellation, the large A contribution to 
the integral is 



Ai(m) 



dA 



,-m\ 2 /2 



(v / 2lrA) 1 - m (-Av / 2) 



2m 4 
A2" ~ A 1 



(65) 
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The first two terms give contributions that are not divergent when m — ¥ 0, hence they are sublcading with respect to 
the 1/A 4 term that gives a finite contribution. We conclude that Ai(m) has a finite limit given by 



A!(0) = 



dA 



0i(A)V _ A 0i(A) 



2-2A 2 -4f?i^V + 6A ^ 



e (A) 



e (A) 



e (A); e (A) 

Instead, the leading large A behavior of the integrand of A2(m) is, changing variable to y = i/m\ 

dA 



(-A\/2) w 1.6 (66) 



A 2 (m) 



2- 



~ mX / 2 (v / 2^A) 1 " m (A 2 - 1) [-m + 



1 



dye- y2/2 y 3 



2 
m 



(67) 



We conclude that y Aq rj 0.8m. Finally, we can show similarly that for small m 

dA 



Gi(m) ~ — m 



e - mA / 2 (\/27A) 1 - m (-Av / 2) 



dye~ v /V 



hence 



1-m 2 V ^ G 



/ 4to 



4 1 



7T 0.8 X M 



1.41m 



-1/2 



Both asymptotic results for Aq and <^g are perfectly consistent with the numerical data. 



E. The 2RSB solution 



(68) 



(69) 



When the 1RSB solution becomes unstable, one must consider further RSB. We performed a 2RSB calculation. 
Then we can linearize the 2RSB solution close to the 1RSB one and obtain the line at which the 2RSB provides a 
better maximization of the entropy, hence becoming stable. This provides an independent calculation of the 1RSB 
instability, which wc verified to be coherent with the one reported above. A complete characterization of the 2RSB 
solution (as well as the 3, 4, ■ ■ • , ooRSB ones) will be presented in future papers of this series. 



VI. COMPUTATION OF THE DYNAMIC EXPONENTS FROM THE CUBIC EXPANSION 

The same strategy allows to obtain the cubic terms in the expansion. From these, following the procedure of 
Ref. [15, 49], one can compute the mean-field dynamical critical exponents at the dynamical glass transition, the 
so-called exponent parameter of Mode-Coupling theory, Amct- Although this calculation is not the main scope of 
this paper, we report it in this section. 

Let us define, following the same notation as for the second order terms (hence for a ^ b, c d, e ^= f which we 
omit from now on) 

W aW = . (70) 

da a<b da c<d oa e< f 

Exploiting the replica symmetry, the two coefficients w\ and W2 can be written in the following form 



W\ — W a b,bc,ca — 3W a b,ac,bd + 3W ac ,bc,de — W a b,cd,ef , 

1 3 „ „„, . „„, „„ r . „„, (71) 

2 W a b,ab,ab ~ 3W a b, a b,ac + ^ 



' 2W a b, c d,ef ) 



and we then have Amct = f^- 
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Defining 



W. 



CO 



5 3 T[i 



ab ^i Sv a<b Sv c<d Sv e<f 



ab ' cd ^ 5a a<b 5a c<d 5a e<f 



0=2a iasB 

logdet(d m ' m ) 



wc have 



W ab , cd , ef = W^l d>ef - 8cpW$ cd!ef , 



hence we have a similar relation for w\ and u>2- We now compute these two terms separately. 



(72) 



(73) 



A. The entropic term 

Following the same strategy as in Sec. VB, we obtain 



W 



(E) 



ab,cd,ef 



Tr 



Sa a 



5a : 



5 a 



e<f 



Tr 



5a a 



Using Eq. (71), the results of Sec. VB and performing the traces we obtain 



5a 



5a 



c<d 



1/' 



(E) 



A 3 



wi E) = 



(74) 



(75) 



B. The interaction term 



The interaction term is 



W, 



5 3 T[i 



ab - cd ^f 5v a<b 8v c<d 5v e<f 



= (f(n a ,n b )f(n c ,n d )f(n e ,nf)) 



(76) 



Using the expression of /, expanding the products, and simplifying many monomials using the symmetries (e.g. 
n a n d n 2 j\ = (n^n^ri c )), we obtain 



(/) _ / n a n^n^. ~ 3 n a n b n cnd + in a n b n c n d n e - n a n b n c n d n e nf 



(I) _ I {l/2)n a nl - in\n\n c + 2n a n b n c n d + (9/2)ngnj;n c n rf - Qn 2 a n b n c n d n e + 2n a n b n c n d n e n f 

n6 



(77) 



Now we convert the average over n into an average over A, using Eq. (48). Performing the derivatives and exploiting 
similar symmetries to simplify the result we obtain 



tr 



(I) 



1 



(2A) 3 
(i) 1/1 



( — 1 + 3Aj — 3A1A2 — 3A^A^ + 6A^A2A3 + AjA^Ag — 3A1A2A3A4 — 3AJA2A3A4 + 3A1A2A3A4A5 — A1A2A3A4A5A6) 

9, 



. — A^Aj — 3AJA2A3 + 2AJA2A3A4 + — A1A2A3A4 — 6A1A2A3A4A5 + 2A1A2A3A4A5A6 ) . 
(2A) 3 \ 2 2 



(78) 
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This result, together with Eq. (53), allows for the explicit computation of these terms. After some simplifications, we 
obtain 

(-0 - 5_/«„m™-7 



< = -7^( e °( A ) I °( A ) + - e (A)9 2 (A)] 2 

(2A) 3 \ 

{(m - 3A 2 )6 (A) 2 + (m - 6)6i(A) 2 + 6 (A)[6A6i(A) - (m - 3)6 2 (A)]} 

w 2 /} = -i-(e (A) m - 7 [2e 1 (A) 3 - 39 (A)e 1 (A)e 2 (A) + e (A) 2 e 3 (A)] ( 7Q ) 

{2A 3 9 (A) 3 + 2(m - 6)6i(A) 3 + 36 (A)e 1 (A)(4Aei(A) - (m - 4)6 2 (A)) 
+ e (A) 2 [-6A(A6i(A) + 6 2 (A)) + (m - 2)6 3 (A)]}) , 

which can be easily computed numerically. 

C. Numerical result 

Collecting all the terms together we obtain the final result 



(/) 



_w 2 _ ~8<M , sm 

AmCT — — — 777 ■ (o(i) 

wi 2/A 3 -8lpw[ T 

When computed at the dynamical transition with m = 1, <p = tpg, = 4.80677 and A = Ad — 0.57668 given by the 
solution of Eq. (25), we obtain 

Amct = 0.70698 , (81) 

which implies that the MCT exponents are a = 0.324016, b = 0.629148 and 7 = 2.33786. The result for 7 is roughly 
consistent with the numerical estimate of Rcf. [41]. 

VII. PHENOMENOLOGICAL EXTENSION TO FINITE DIMENSIONS 

We can obtain quantitative results in finite d by a phenomenological extension of Eq. (15). First we go back to 
non-rescaled density and we rearrange it as 

s[a] = 1 - log p - 2 d ~V + dlogm + (m ~ ^ \og{2neD 2 / d 2 ) + | logdct(a m ' m ) + 2 d ~ V [1 -^(26)] . (82) 

We now recognize that su q = 1 — logp — 2 d ~ 1 (p. Furthermore, by comparison with the finite d results obtained in 
the small cage expansion [9] , we know that the interaction term is renormalized by the contact value of the liquid 
correlation yu q (tp). We therefore can propose the following form for the entropy: 

s[&] = s Uq (<p) + dlogm + (m 2 1)d log(27reD 2 /d 2 ) + I logdet(d m < m ) + 2 d " Vj/« g (</?) [1 - F {2a)] . (83) 
In the 1RSB scheme we obtain 

rTi , , d, (m-l)d (m-l)d, (2itD 2 A\ h _, , . „ 
s[A]= Sltq { V ) + -\ogm+± + ~ log I d2 1 +2 d ^ViMQmiA) , 

F m {A) , (84) 



A K = -4-2%M Am(4 ) 

Although these equations are not obtained from a consistent derivation, recalling that for small A we have G m {A) ~ 
\/C4Gi(m) and Gi(m) = 2<5o(w), they reproduce the small cage expansion of Ref. [9] at the leading order in A. Note 



21 



d 


<y5G 


<PGCP 


PG 


3 


0.683581 


0.683657 


26727 


4 


0.486755 


0.486874 


16374 


5 


0.330586 


0.330718 


12535 


6 


0.218074 


0.218203 


10119 


7 


0.140074 


0.140189 


8469 


<s 


0.0876190 


0.0877137 


7407 


9 


0.0534490 


0.0535198 


6805 


10 


0.0318889 


0.0319377 


6550 


11 


0.0186760 


0.0187075 


6548 


12 


0.0107756 


0.0107949 


6724 


13 0.00614419 


0.00615559 


7019 



TABLE I: Values of f>GCP [9] and of <^G and pc for several dimensions. Note that these values correspond to the equilibrium 
Gardner transition. For out-of-equilibrium states (which are the one produced in all experiments and numerical simulations), 
we expect that the Gardner instability will happen at lower, possibly much lower, pressures and densities. 



that when expressed in terms of A and m, the equation for the stability A# = is exactly the same as in d — > oo, 
Eq. (58). Hence the result for Aq(itl) is independent of dimension. 

We will check a posteriori that even in d = 3 the Gardner transition happens at very large pressure, hence m and A 
are small. So we can use the asymptotic expansions to obtain quantitative estimates. The procedure is the following: 

• Recall that at small m wc have y Aq 0.8m. 

• Now we obtain tpQ{rn) (or better mc(p)) by solving 



Gi(m) 



1 — m 



0.4m w 0.71-v/m 



We recall from the analysis of Rcf. [9] that 
m* ~ n ((fccp - <p) , 



m G 



0.71 x 2 d <py liq ((p) 



id—l 



yuq{y) - d 



Vliq(<P) 



1-d 

9 



V=VGCP 



(85) 



(86) 



• The Gardner transition happens when the two lines cross, hence (pc is the solution of mcfip) = m*(ip) which 
can be easily found numerically once an equation of state for the liquid has been chosen. Here we use the 
Carnahan-Starling equation already used in Ref. [9] . 



Finally we use the result [9] 



p(<p, 



G) 



dlfiGCP 
fGCP - V 



(87) 



to estimate the Gardner pressure pa = PifG)- 



The numerical values of the Gardner pressure are reported in Tab. I. Note that the non-monotonicity of pc at low 
dimension could be an artifact of the approximations used above. Note also that pc is always much larger than the 
pressure at the glass transition (reported in Rcf. [9, 43]). Still, the reader should keep in mind that the value of pc 
corresponds to the Gardner instability of the equilibrium ("ideal") glass. According to the phase diagram of Fig. 1, 
we expect that the Gardner instability for the metastable states that are reached out of equilibrium will happen at 
lower pressures. Unfortunately, quantifying this effect requires the use of "state following" techniques [50] and goes 
far beyond the scope of this article. 

The limit d — > oo is recovered as follows. Recall that (pgcp ~ logd [9]. Moreover, yu q — > 1. Hence [i ~ 2 d ~ 1 /d and 
the equation for (pG becomes 



\($gcp - <Pg) = ( n =. z ~ 
2 \0.71 x IpG 



1 



0.71 x ipgcp 



1 



0.71 x logd 



(88) 
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which shows that the distance between <pc and <pgcp shrinks as (logd) -2 and the Gardner pressure diverges as 
Pg ~ d(logcf) 3 , as was sketched in Fig. 1. 

Finally, at this level of approximation, it can be easily shown that Amct does not depend on dimension. The small 
dependence of Amct reported in Ref. [41] should be explained by corrections to this approximation, that have been 
neglected here. 

VIII. CONCLUSIONS 

In this paper we were able to investigate the possibility of a Gardner transition for hard spheres in large spatial 
dimensions. Such a study was never been done before for particle systems, and was possible thanks to the expression 
of the entropy in terms of the overlap matrix obtained in the first paper of this series, and extended here to obtain 
Eqs. (15) and (16): because this expression has been shown to be exact, any discrepancy or instability is only 
attributable to the instability of the 1RSB ansatz. 

The 1RSB solution is unstable in the equilibrium glass phase for (reduced) pressure higher than the Garner pressure 
Pg [12], and for metastable states in large region of the phase diagram, just as in p-spin glasses [21, 22]. We provided 
an estimate of the Gardner pressure in finite dimensions, finding that it is quite high; moreover we showed that pg 
diverges (slowly) with increasing dimension. These pressures are directly accessible to numerical simulations, hence we 
expect that this transition should be quite easy to detect numerically. Estimating analytically the transition pressure 
for metastable glasses would be very useful to guide numerical simulations: this is however hard, as it requires a state 
following computation [50] . Although this is possible in principle, we leave it for future work. We expect that in any 
case the instability will happen at lower pressure for metastable glasses than for the ideal glass. 

The physical consequences of this instability are very intriguing but for the moment not all its implications have 
been worked out. In fact, even for the simplest p-spin glasses, the impact of the Gardner instability on the out- 
of-equilibrium dynamics is not completely understood from a technical point of view [22]. The structure of the 
metastable states of the p-spin glass model and its impact on the out-of-equilibrium dynamics are being actively 
investigated [12, 21, 22, 25, 50-53] and making progress on this simpler model will be crucial for understanding 
the technically more involved hard sphere case. We expect (hope) that the scenario we proposed in Sec. II will be 
confirmed by these studies. 

Let us recall here some speculations on the possible impact of the Gardner instability on the physics of jamming 
that we discussed in this paper, leaving a more detailed investigation for future works. 

• It is reasonable to expect that at the Gardner transition, the 1RSB solution will transform continuously into 
a full RSB solution, although we have not yet constructed this explicitly. Such a solution describes a situation 
where glassy states are arranged in a complex and correlated pattern [6] . More importantly, they are marginally 
stable [6, 54]. This means that the spectrum of vibrations around a glassy state displays many soft modes. 
Hence, it is likely that a full RSB description of the problem will allow one to obtain information on the soft 
modes that are observed at the jamming transition [13, 14, 55], especially those of frequency below the gap lo* . 
Some steps in this direction have been already performed in Ref. [40], where the analogy with the full RSB 
physics of the Shcrrington-Kirkpatrick model was noticed. 

• There should be several signatures of the Gardner transition. Suppose that the hard sphere system is prepared 
in a glass state in the region where the 1RSB solution is stable, and that pressure in slowly increased approaching 
the instability. As mentioned in Sec. II, the spin glass susceptibility [6] (which in this context is a four-point static 
susceptibility) diverges on approaching the instability. Moreover, even if the system was already equilibrated in 
the initial glass state, aging effects should appear below the instability when the state breaks down into many 
correlated sub-states. 

• The aging curves, and in particular the fluctuation-dissipation plots, should give a good indication of the 
transition: at pressures above the Gardner pressure these plots should crossover from two straight lines to two 
straight lines joined by a curved segment (although detecting the curved part might be numerically challenging). 

• There is probably a relation between the Gardner transition and the "dynamic criticality" defined in Ref. [33]. 
As noticed there, all the anomalous scalings at the jamming transition are related to the scaling of the cage 
radius with pressure, A ~ p~ 3 / 2 . Hopefully, this scaling, which is not found in the 1RSB solution, could be a 
property of the full RSB phase. In fact, it is well known that in the Sherrington-Kirkpatrick the presence of a 
full RSB phase changes the scaling of the overlap at low temperatures. 

• Brito and Wyart [27] have demonstrated that close to jamming, the dynamics is characterized by sudden 
"cracks" at which the system leaves abruptly a locally stable structure to find a new one. At these cracks, 
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the displacements of the particles are strongly correlated with the lowest frequency eigenvectors of the stability 
matrix of the structure that the system is leaving. As mentioned in the introduction, this is not what is expected 
in a 1RSB phase, where states are locally stable and one has to cross a barrier to jump from one state to the 
other - hence the vibrations at the bottom of the well should give no information on the shape of the barrier. 
However, in a full RSB phase the dynamics is much different and similar to the one found in Ref. [27]. It would 
be very nice to check whether the results of Brito and Wyart really fit into a full RSB picture, for example by 
calculating the spatial distribution of the displacements between two nearby states. 

• Finally, the response of full RSB magnetic systems to an external perturbation is very complex, being character- 
ized by avalanches and intermittency, see e.g. Ref. [56]. This is due to the existence of (relatively low) barriers 
separating nearby states - all this within a large basin. By analogy, we would expect the response of a hard 
sphere system to a mechanical perturbation in the full RSB phase to be similarly complex. Hence, the Theolog- 
ical properties in this phase could be very interesting and could explain some of the anomalous behavior found 
around the jamming transition. Analytical computations might be possible following the strategy introduced in 
Ref. [57-59]. 

From the technical point of view, the next step to make progress is to investigate the i^RSB solutions, with K = 
2, 3, 4, • • • , eventually with K — > oo that corresponds to full RSB. Following Gardner's example in the p-spin model, 
this can be done just below the Gardner transition. This investigation is in progress and will provide some answers 
to the above questions. In parallel, numerical simulations should be performed to detect the 1RSB instability. Also, 
an exact solution of the dynamics, along the lines of Ref. [60], could provide very useful complementary informations. 

To conclude, let us mention that in this paper, from the expansion of the cubic terms around the 1RSB solution 
(see Section VI), we obtained an estimate of the mean- field dynamical critical exponents at the dynamical transition 
(the so-called Mode-Coupling Theory exponent parameter Amct) [61]. We found that for d — > oo, Amct = 0.70698, 
which is consistent with numerical simulations [41]. This is important because a previously attempted calculation 
from the replicated HNC equations [15] gives results that are quantitatively bad. The fact that in d — > oo we can 
obtain a good result implies that the negative result obtained in Ref. [15] has to be attributed to the poor quantitative 
performances of the replicated HNC approximation, which indeed were already known [9] . Unfortunately, also the 
approach presented here gives poor quantitative results for the dynamical glass transition in low dimensions [9]. 
Obtaining an accurate theory of the dynamical glass transition in low dimension by improving the replicated HNC is 
therefore very important, see [62] for a preliminary step in this direction. 
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Appendix A: Computation of the Jacobian J(q) 
The Jacobian J(q) is defined in the following way 



IP 



.6=1 



d d u 1 . . .d d « m -i S(q ab 



(Al) 



,<b 



Let us introduce the following notation. We define adx (m — 1) matrix U whose first column contains the components 
of the d-dimensional vector u\ , the second column contain the components of 1*2 and the m — 1 column contains the 
components of u m -i- Moreover we define the m x m matrix q and its reduced version that is the matrix q( m > m ) which 
is obtained from q by deleting the last column and the last row. Let us consider the integral in the expression (Al): 



d d u 1 ...d d u m - 1 8(q ab 



i<b 



dUS 



(m,m) 



(A2) 
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where the last integral is with a flat measure over the entries of the matrix U. To perform this computation we start 
from a simple case. Let us consider the case where the matrix q( m > m ) has a diagonal structure 



q {m ' m) = to 



(to)afc — QaaSab 



This means that 



dUS 



m—1 



(l,m-l 



<7 (m ' m) -U T U = dV . . . d d u m ^ J] 5 (q aa - \u a \ 2 ) S(u a ■ u b ) 



(A3) 



(A4) 



a=l 



a<b 



The second Dirac delta function says that the vectors u have to be all orthogonal one to the other. The first one tells 
us about the length of these vectors. By going to polar coordinates we have that the previous expression is given by 



d d Ul . . . d d 



l^m— 1 



d-1 



m—1 



du m _iu d l } 1 Y[ S (q aa - u 2 a ) 5{u a u b u a ■ u b ) 



(A5) 



where we have denoted with d d u the angular integration in d dimensions and with u a the unit vector parallel to u a . 
Rewriting the Dirac deltas in different ways the previous expression can be written in the following way 



/r°° t°° / ra 1 i 

d% . . .d d u m ^ / d^" 1 . . . / du m ^u d -\ ] I -—5 (, 
JO JO \o=l V^ aa 



Wo) X 



n 

a<b 



1 



^Jqaa^bb 



1,-m-l 

5{iL a ■ Ub) 



a<b 



where we have 



C m .d — 2 



a=l 



l.m-1 



deta (m ' m) 



(d-m)/2 



(A6) 



d d ui . . . d d w m _i 



n *(< 



(A7) 



a<b 



The integral that appears in the expression for C m> d is very simple. In fact the Dirac deltas tell us that the unit 
vectors u must be all orthogonal. So we have to compute the phase space accessible to them. This can be done 
iteratively. Suppose that we have only one unit free vector. It has a phase space available given by the solid angle in 
d dimensions which is fid- Then we add a second unit vector orthogonal to the first one. Clearly it has a phase space 
available that is the solid angle in the space orthogonal to the first vector that is Qd-i- Going on by iteration we have 



C m ,d — 2 m OdQd_i . . . Q.d-m+2- 



(A8) 



Now we want to generalize this to a matrix g( m > m ) that is not diagonal. Because the matrix g( m ' m ) is symmetric we 
can diagonalize it and we can write 



~(m,m) = A -iq oA ) dct A = 1 , dct q D = det q^^ , A T = A" 1 



It follows that 



dUS 



_ jjTjj 



If we change integration variables 



= / dUS [A -1 $ D A - U T U] = / dUS [q D - AJ7 T [/A _1 ] 



U' = UA~ 



(A9) 
(A10) 
(All) 



that has a Jacobian which is unitary due to the orthogonality property of the matrix A, then the previous expression 
becomes 



J dUS [q D - A[/ T [/A _1 ] = J dUS [q D - U T U] , 



that is the same that we discussed above. It follows that 



J(q) = m d C m4 J] S 



a=l Lb=l 



det g (m ' m) 



(d-m)/1 



(A12) 



(A13) 



which is the result that was used in Ref. [11]. 
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